####
## Project:            Do Corporate Regulations Deter or Stimulate Investment? The Effect of the OECD Anti-Bribery Convention on FDI
## File name:          02-firm_analysis.R
## Description:        Analyses firm-level data and exports tables and plots with results
## Programmed by:      Lorenzo Crippa
####

# Setup ----
graphics.off()
rm(list = ls())
setwd("") # INSERT THE REPLICATION FOLDER FILE PATH HERE 
set.seed(1234)

library(tictoc)
library(arm) 
library(lme4)
library(showtext)
library(modelsummary)
library(countrycode)
library(tidyverse)
library(tibble)
source("code/aux_functions.R")

# plots' font:
font_add("Times New Roman", "times new roman.ttf")
showtext_auto()

# simplify modelsummary LaTeX format:
options("modelsummary_format_numeric_latex" = "plain")

theme <- theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
               panel.background = element_rect(fill = "transparent", color = NA),
               plot.background = element_rect(fill = "transparent", color = NA),
               axis.line = element_line(colour = "black", linewidth = 0.3),
               text = element_text(colour = "black", size = 14, family = "Times New Roman"))
theme_set(theme)

# Import data ----
orbis <- read_rds("data/firm_level.rds")

# Drop observations not needed ----
# drop observations where isocodeHome == isocodeHost and those where hosttot == 0
orbis <- orbis %>% filter(isocodeHome != isocodeHost)
orbis <- orbis %>% filter(hosttot > 0)

# Summary statistics ----
firm_level <- orbis %>%
  group_by(parent) %>%
  select(isocodeHome, oecdABC, lgdpHome, gdpgrowHome, ljiHome,
         lage, lassets, lhosts, naics3, naics2, nace, 
         industry_bribery, home_bribery) %>%
  distinct()

## Table C.1 ----
datasummary_balance(formula = ~oecdABC,
                    data = firm_level %>%
                      ungroup() %>%
                      select(-parent, -isocodeHome, -naics3, -naics2, -nace) %>%
                      mutate(oecdABC = factor(oecdABC, levels = c(0, 1),
                                              labels = c('Control', 'Treatment'))) %>%
                      relocate(lage, lassets, lhosts, industry_bribery) %>%
                      rename(`Home GDP (log)` = "lgdpHome",
                             `Home GDP Growth (%)` = "gdpgrowHome",
                             `Home Judiciary Indep.` = 'ljiHome',
                             `Firm Age (log)` = 'lage',
                             `Firm Assets (log)` = 'lassets',
                             `Firm Host Countries (log)` = 'lhosts',
                             `Industry with enforcement` = 'industry_bribery'),
                    fmt = 2, dinm_statistic = "std.error",
                    align = 'lcccccc',
                    output = 'tables/supp_tab_c1.tex',
                    title = 'Description of covariates for treated and control firms')

## Figure B.1 ----
firm_level %>%
  mutate(naics2 = as_factor(naics2)) %>%
  group_by(naics2, oecdABC) %>% 
  summarise(count_occurr = n()) %>%
  group_by(oecdABC) %>%
  mutate(tot = sum(count_occurr),
         perc = 100*count_occurr/tot,
         oecdABC = factor(oecdABC, levels = c(0, 1),
                          labels = c('Control', 'Treatment'))) %>%
  ggplot(aes(x = reorder(naics2, -perc), y = perc)) + 
  geom_col(fill = "grey60") + 
  coord_flip() + ylab("%") + xlab("") + 
  facet_wrap('oecdABC')
ggsave("plots/supp_fig_c1.pdf", device = "pdf", width = 7, height = 5)

## Figure B.2 ----
firm_level %>%
  mutate(isocodeHome = as_factor(isocodeHome)) %>%
  group_by(isocodeHome, oecdABC) %>% 
  summarise(count_occurr = n()) %>%
  group_by(oecdABC) %>%
  mutate(tot = sum(count_occurr),
         perc = 100*count_occurr/tot,
         oecdABC = factor(oecdABC, levels = c(0, 1),
                          labels = c('Control', 'Treatment')),
         isocodeHome = countrycode(isocodeHome, origin = 'iso2c', destination = 'country.name.en')) %>%
  ggplot(aes(x = reorder(isocodeHome, -perc), y = perc)) + 
  geom_col(fill = "grey60") + 
  coord_flip() + ylab("%") + xlab("") + 
  facet_wrap('oecdABC', scales = 'free_y')
ggsave("plots/supp_fig_c2.pdf", device = "pdf", width = 7, height = 6)

## Table C.2 ----
# variables not to show in the table
noshow <- c("dyad", "isocodeHome", "isocodeHost",
            "naics3","naics2", "nace", 
            "hosttot", 
            "gdppcHome", "fuelgdp",
            "oecdHome", "oecdHost"
)

dat_summary <- orbis %>%
  select(-all_of(noshow)) %>%
  rename(`Subsidiary` = "subsid",
         `Subsidiary (extended)` = "subsidExtend",
         `OECD Ratifier` = "oecdABC",
         `Host PACI` = "paciHost2005",
         `Host PACI (2012)` = "paciHost",
         `Host Bribery V-Dem` = "v2exbribeHost",
         `Host CCE` = "cceHost",
         `Industry with enforcement` = 'industry_bribery',
         `Host with enforcement` = 'host_bribery',
         `Home with enforcement` = 'home_bribery',
         `Home GDP (log)` = "lgdpHome",
         `Home GDP Growth (%)` = "gdpgrowHome",
         `Home Judiciary Indep.` = 'ljiHome',
         `Host GDP (log)` = "lgdpHost",
         `Host FDI (GDP %)` = "fdigdpHost",
         `Host GDP per capita` = "gdppcHost",
         `Host GDP Growth (%)` = "gdpgrowHome",
         `Host Judiciary Indep.` = 'ljiHost',
         `Host Trade (GDP %)` = "tradeHost",
         `Host Democracy` = "demHost",
         `Host POLCON III` = "pcHost",
         `Dyad Distance (hundreds of km)` = "dist",
         `Dyad Common Language` = "comlang",
         `Dyad Colonial Relation` = "colony",
         `Dyad BIT` = 'bit',
         `Firm Age (log)` = 'lage',
         `Firm Assets (log)` = 'lassets',
         `Firm Host Countries (log)` = 'lhosts')

datasummary(data = dat_summary,
            formula = All(dat_summary) ~ N + Mean + SD + Min + Max,
            fmt = 3,
            output = "tables/supp_tab_c2.tex",
            title = "Firm-level data. Summary statistics")

# Recenter covariates for models ----
# recentering continuous control variables at means to help convergence
orbis$ljiHome 		<- recenter(orbis$ljiHome)
orbis$ljiHost 		<- recenter(orbis$ljiHost)
orbis$pcHost 	<- recenter(orbis$pcHost)
orbis$demHost	<- recenter(orbis$demHost)
orbis$lassets		<- recenter(orbis$lassets)
orbis$lage		<- recenter(orbis$lage)
orbis$gdpgrowHome	<- recenter(orbis$gdpgrowHome)
orbis$lgdpHost 	<- recenter(orbis$lgdpHost)
orbis$lgdpHome 	<- recenter(orbis$lgdpHome)
orbis$gdppcHost 	<- recenter(orbis$gdppcHost)
orbis$tradeHost	<- recenter(orbis$tradeHost)
orbis$fdigdpHost	<- recenter(orbis$fdigdpHost)
orbis$lhosts	<- recenter(orbis$lhosts)
orbis$dist		<- recenter(orbis$dist)
orbis$fuelgdp <- recenter(orbis$fuelgdp)
orbis$gdppcHome <- recenter(orbis$gdppcHome)

# Main text analysis ----
tic()
base <- glmer(subsid ~ oecdABC * paciHost2005 +
                oecdABC * I(paciHost2005*paciHost2005) +
                (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
              data = orbis,
              family = binomial(link="logit"), nAGQ = 0)
summary(base)
toc() # 17 seconds

tic()
country <- glmer(subsid ~ oecdABC * paciHost2005 +
                   oecdABC * I(paciHost2005*paciHost2005) +
                   lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                   ljiHost + pcHost + demHost +
                   lgdpHome + gdpgrowHome +
                   ljiHome + 
                   (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
                 data = orbis,
                 family = binomial(link="logit"), nAGQ = 0)
summary(country)
toc() # 17 seconds

tic()
dyad <- glmer(subsid ~ oecdABC * paciHost2005 +
                oecdABC * I(paciHost2005*paciHost2005) +
                lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                ljiHost + pcHost + demHost +
                lgdpHome + gdpgrowHome +
                ljiHome + 
                bit + comlang + colony + dist +
                (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
              data = orbis,
              family = binomial(link="logit"), nAGQ = 0)
summary(dyad)
toc() # 29 seconds

tic()
full <- glmer(subsid ~ oecdABC * paciHost2005 +
                oecdABC * I(paciHost2005*paciHost2005) +
                lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                ljiHost + pcHost + demHost +
                lgdpHome + gdpgrowHome +
                ljiHome + 
                bit + comlang + colony + dist +
                lassets + lage + lhosts +
                (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
              data = orbis,
              family = binomial(link="logit"), nAGQ = 0)
summary(full)
toc() # 33 seconds

tic()
full.ind <- glmer(subsid ~ oecdABC * paciHost2005 +
                    oecdABC * I(paciHost2005*paciHost2005) +
                    lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                    ljiHost + pcHost + demHost +
                    lgdpHome + gdpgrowHome +
                    ljiHome + 
                    bit + comlang + colony + dist +
                    lassets + lage + lhosts +
                    (1|isocodeHome) + (1|isocodeHost) + (1|naics2)  + (1|dyad),
                  data = orbis,
                  family = binomial(link="logit"), nAGQ = 0)
summary(full.ind)
toc() # 54 seconds

tic()
full.ind.ex <- glmer(subsidExtend ~ oecdABC * paciHost2005 +
                       oecdABC * I(paciHost2005*paciHost2005) +
                       lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                       ljiHost + pcHost + demHost +
                       lgdpHome + gdpgrowHome +
                       ljiHome + 
                       bit + comlang + colony + dist +
                       lassets + lage + lhosts +
                       (1|isocodeHome) + (1|isocodeHost) + (1|naics2)  + (1|dyad),
                     data = orbis,
                     family = binomial(link="logit"), nAGQ = 0)
summary(full.ind.ex)
toc() # 42 seconds

## Table 1 ----
# note: some minor manual aesthetic tweaks to the table columns in the paper and appendix
cf_map <- c("oecdABCI(paciHost2005 * paciHost2005)" = "OECD Ratifier X Host PACI sq.", 
            "oecdABCpaciHost2005" = "OECD Ratifier X Host PACI",
            "oecdABC" = "OECD Ratifier",
            "I(paciHost2005 * paciHost2005)" = "Host PACI sq.",
            "paciHost2005" = "Host PACI")

rows <- tribble(~term, ~base, ~country, ~dyad, ~full, ~full.ind, ~full.ind.ex,
                "Random intercepts", "Yes", "Yes", "Yes", "Yes", "Yes", "Yes",
                "Country-level contorls", "", "Yes", "Yes", "Yes", "Yes", "Yes",
                "Dyad-level controls", "", "", "Yes", "Yes", "Yes", "Yes",
                "Firm-level controls", "", "", "", "Yes", "Yes", "Yes",
                "Industry intercepts", "", "", "", "", "Yes", "Yes",
                "Extended dep. var. data", "", "", "", "", "", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(base@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(country@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(dyad@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.ind@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.ind.ex@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(base@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(country@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(dyad@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.ind@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.ind.ex@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(base)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(country)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(dyad)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.ind)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.ind.ex)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(11:18, 20)

modelsummary(list("(1)" = base,
                  "(2)" = country,
                  "(3)" = dyad,
                  "(4)" = full,
                  "(5)" = full.ind,
                  "(6)" = full.ind.ex),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/main_tab_1.tex",
             title = "Firm-level data. The effect of the OECD Convention on probability of subsidiary incorporation. Multilevel logit models"
)

## Table B.3 ----
cf_map <- c("oecdABCI(paciHost2005 * paciHost2005)" = "OECD Ratifier X Host PACI sq.", 
            "oecdABCpaciHost2005" = "OECD Ratifier X Host PACI",
            "oecdABC" = "OECD Ratifier",
            "I(paciHost2005 * paciHost2005)" = "Host PACI sq.",
            "paciHost2005" = "Host PACI",
            "lgdpHost" = "Host GDP (log)", 
            "gdppcHost" = "Host GDP per capita",
            "fdigdpHost" = "Host FDI (GDP perc.)", 
            "tradeHost" = "Host Trade (GDP perc.)", 
            "ljiHost" = "Host Judiciary Indep.", 
            "pcHost" = "Host POLCON III", 
            "demHost" = "Host Democracy",
            "lgdpHome" = "Home GDP (log)", 
            "gdpgrowHome" = "Home GDP Growth (perc.)", 
            "ljiHome" = "Home Judiciary Indep.",
            "bit" = "Dyad BIT", 
            "comlang" = "Dyad Common Language", 
            "colony" = "Dyad Colonial Relation", 
            "dist" = "Dyad Distance",
            "lassets" = "Firm Assets (log)", 
            "lage" = "Firm Age (log)", 
            "lhosts" = "Firm Host Countries (log)",
            "(Intercept)" = "Constant")

rows <- tribble(~term, ~base, ~country, ~dyad, ~full, ~full.ind,
                "Random intercepts", "Yes", "Yes", "Yes", "Yes", "Yes",
                "Industry intercepts", "", "", "", "", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(base@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(country@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(dyad@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.ind@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(base@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(country@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(dyad@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.ind@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(base)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(country)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(dyad)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.ind)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(47:50, 52)

modelsummary(list("(1)" = base,
                  "(2)" = country,
                  "(3)" = dyad,
                  "(4)" = full,
                  "(5)" = full.ind),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/supp_tab_d1.tex",
             title = "Firm-level data. The effect of the OECD Convention on probability of subsidiary incorporation. Multilevel logit models (full disclosure)"
)

## Figure 3 ----
# generate the plot and divide in three panels to help the discussion
ci <- .95
r <- seq(min(orbis$paciHost2005, na.rm = T), max(orbis$paciHost2005, na.rm = T), by = .025)

# calculate effect:
change <- marg_eff(base, range = r, ci.lvl = ci)

# add labels for number of hosts:
n_a <- length(unique(orbis$isocodeHost[orbis$paciHost2005 <= 3.3 & is.na(orbis$paciHost2005) == FALSE]))
n_b <- length(unique(orbis$isocodeHost[orbis$paciHost2005 > 3.3 & orbis$paciHost2005 <= 5.98 & is.na(orbis$paciHost2005) == FALSE]))
n_c <- length(unique(orbis$isocodeHost[orbis$paciHost2005 > 5.98 & is.na(orbis$paciHost2005) == FALSE]))

# because I want to plot result with distribution of moderator at the bottom of the plot, I need to shift up all the effects on the Y
shift <- 100 # this is the factor by which to shift estimates up, so that we can place the histogram at the bottom
enlarger <- 200 # this is the factor to multiply histogram and make it larger (or it won't be visible...)

change %>% 
  ggplot(aes(x = x, y = change+shift)) + geom_line() +
  geom_ribbon(aes(ymin = lower+shift, ymax = upper+shift, alpha = I(.1))) +
  geom_histogram(data = orbis %>%
                   select(oecdABC, paciHost2005, dyad, isocodeHost, isocodeHome) %>%
                   filter(complete.cases(.)),
                 aes(x = paciHost2005, y = enlarger*((after_stat(count))/sum(after_stat(count))))) +
  geom_hline(yintercept = 0+shift, color = I("gray40"), size = 0.3) +
  geom_vline(xintercept = 3.3, color = I("gray40"), linetype = I("dashed"), linewidth = 0.3) +
  geom_vline(xintercept = 5.98, color = I("gray40"), linetype = I("dashed"), linewidth = 0.3) +
  geom_label(aes(label = "(a)", y = 80+shift, x = 1.5), family = "Times New Roman", colour = "black", 
             label.size = NA) +
  geom_label(aes(label = "(b)", y = 80+shift, x = 4.5), family = "Times New Roman", colour = "black", 
             label.size = NA) +
  geom_label(aes(label = "(c)", y = 80+shift, x = 7.5), family = "Times New Roman", colour = "black", 
             label.size = NA) +
  geom_label(aes(label = paste0("N = ", n_a), 
                 y = 70+shift, x = 1.5), family = "Times New Roman", colour = "black", 
             label.size = NA) +
  geom_label(aes(label = paste0("N = ", n_b), 
                 y = 70+shift, x = 4.5), family = "Times New Roman", colour = "black", 
             label.size = NA) +
  geom_label(aes(label = paste0("N = ", n_c), 
                 y = 70+shift, x = 7.5), family = "Times New Roman", colour = "black", 
             label.size = NA) +
  scale_y_continuous(breaks = c(shift-50, shift, shift+50), labels = function(x) paste0(x-shift, "%")) + 
  xlab("Host PACI") + ylab("Change in predicted probability of subsidiary incorporation")
ggsave("plots/main_fig_3.pdf", device = 'pdf', width = 5, height = 6)

max(change$change) # +27%
min(change$change) # -52%

# what is the probability of investment for unregulated companies in countries at the peak of the empowerment effect?
mean(orbis$subsid[orbis$paciHost2005 <= 3.2 & orbis$paciHost2005 >= 2.8 & !is.na(orbis$paciHost2005)], na.rm = TRUE) # 0.02

# what is the probability of investment for unregulated companies in countries at the peak of the deterrence effect?
mean(orbis$subsid[orbis$paciHost2005 >= 7.5 & !is.na(orbis$paciHost2005)], na.rm = TRUE) # 0.001

# Robustness tests ----
## Corruption measures/specifications ----
tic()
full.order1 <- glmer(subsid ~ oecdABC * paciHost2005 +
                       lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                       ljiHost + pcHost + demHost +
                       lgdpHome + gdpgrowHome +
                       ljiHome + 
                       bit + comlang + colony + dist +
                       lassets + lage + lhosts +
                       (1|isocodeHome) + (1|isocodeHost) + (1|naics2)  + (1|dyad),
                     data = orbis,
                     family = binomial(link="logit"), nAGQ = 0)
summary(full.order1)
toc() # 47 seconds

tic()
full.vdem <- glmer(subsid ~ oecdABC * v2exbribeHost +
                     I(v2exbribeHost*v2exbribeHost)*oecdABC +
                     lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                     ljiHost + pcHost + demHost +
                     lgdpHome + gdpgrowHome +
                     ljiHome + 
                     bit + comlang + colony + dist +
                     lassets + lage + lhosts +
                     (1|isocodeHome) + (1|isocodeHost) + (1|naics2)  + (1|dyad),
                   data = orbis,
                   family = binomial(link="logit"), nAGQ = 0)
summary(full.vdem)
toc() # 59 seconds

tic()
full.CCE <- glmer(subsid ~ oecdABC * cceHost +
                    I(cceHost*cceHost)*oecdABC +
                    lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                    ljiHost + pcHost + demHost +
                    lgdpHome + gdpgrowHome +
                    ljiHome + 
                    bit + comlang + colony + dist +
                    lassets + lage + lhosts +
                    (1|isocodeHome) + (1|isocodeHost) + (1|naics2)  + (1|dyad),
                  data = orbis,
                  family = binomial(link = "logit"), nAGQ = 0)
summary(full.CCE)
toc() # 80 seconds

tic()
full.PACI2 <- glmer(subsid ~ oecdABC * paciHost +
                      I(paciHost*paciHost)*oecdABC +
                      lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                      ljiHost + pcHost + demHost +
                      lgdpHome + gdpgrowHome +
                      ljiHome + 
                      bit + comlang + colony + dist +
                      lassets + lage + lhosts +
                      (1|isocodeHome) + (1|isocodeHost) + (1|naics2)  + (1|dyad),
                    data = orbis,
                    family = binomial(link="logit"), nAGQ = 0)
summary(full.PACI2)
toc() # 53 seconds

### Table E.1 ----
cf_map <- c("oecdABC" = "OECD Ratifier",
            "oecdABCpaciHost2005" = "OECD Ratifier X Host PACI",
            "oecdABCv2exbribeHost" = "OECD Ratifier X Host Bribery V-Dem",
            "oecdABCI(v2exbribeHost * v2exbribeHost)" = "OECD Ratifier X Host Bribery V-Dem sq.",
            "oecdABCcceHost" = "OECD Ratifier X Host CCE",
            "oecdABCI(cceHost * cceHost)" = "OECD Ratifier X Host CCE sq.",
            "oecdABCpaciHost" = "OECD Ratifier X Host PACI (2012)",
            "oecdABCI(paciHost * paciHost)" = "OECD Ratifier X Host PACI sq. (2012)",
            "paciHost2005" = "Host PACI",
            "v2exbribeHost" = "Host Bribery V-Dem",
            "I(v2exbribeHost * v2exbribeHost)" = "Host Bribery V-Dem sq.",
            "cceHost" = "Host CCE",
            "I(cceHost * cceHost)" = "Host CCE sq.",
            "paciHost" = "Host PACI (2012)",
            "I(paciHost * paciHost)" = "Host PACI sq. (2012)",
            "lgdpHost" = "Host GDP (log)", 
            "gdppcHost" = "Host GDP per capita",
            "fdigdpHost" = "Host FDI (GDP perc.)", 
            "tradeHost" = "Host Trade (GDP perc.)", 
            "ljiHost" = "Host Judiciary Indep.", 
            "pcHost" = "Host POLCON III", 
            "demHost" = "Host Democracy",
            "lgdpHome" = "Home GDP (log)", 
            "gdpgrowHome" = "Home GDP Growth (perc.)", 
            "ljiHome" = "Home Judiciary Indep.",
            "bit" = "Dyad BIT", 
            "comlang" = "Dyad Common Language", 
            "colony" = "Dyad Colonial Relation", 
            "dist" = "Dyad Distance",
            "lassets" = "Firm Assets (log)", 
            "lage" = "Firm Age (log)", 
            "lhosts" = "Firm Host Countries (log)",
            "(Intercept)" = "Constant")

rows <- tribble(~term, ~full.order1, ~full.vdem, ~full.CCE, ~full.PACI2,
                "Random intercepts", "Yes", "Yes", "Yes", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(full.order1@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.vdem@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.CCE@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.PACI2@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(full.order1@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.vdem@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.CCE@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.PACI2@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(full.order1)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.vdem)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.CCE)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.PACI2)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(67:69, 71)

modelsummary(list("(1)" = full.order1,
                  "(2)" = full.vdem,
                  "(3)" = full.CCE,
                  "(4)" = full.PACI2),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/supp_tab_e1.tex",
             title = "Firm-level data. Robustness tests of multilevel logit models to different corruption specifications")

# correlation between Host PACI 2005 and Host PACI 2012? (Needed in the text of the appendix)
corr <- cor(orbis$paciHost, orbis$paciHost2005, use = 'pairwise.complete.obs')
corr # 0.98

## Further specification and exclusion of data ----
tic()
original.BB <- glmer(subsid ~ oecdABC * paciHost2005 +
                       oecdABC * I(paciHost2005*paciHost2005) +
                       ljiHost * ljiHome +
                       lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                       ljiHost + pcHost + demHost +
                       lgdpHome + gdpgrowHome +
                       ljiHome + 
                       bit + comlang + colony + dist +
                       lassets + lage + lhosts +
                       (1|isocodeHome) + (1|isocodeHost) + (1|naics2) + (1|dyad),
                     data = orbis,
                     family = binomial(link="logit"), nAGQ = 0)
summary(original.BB)
toc() # 52 seconds

tic()
no.cn <- glmer(subsid ~ oecdABC * paciHost2005 +
                 I(paciHost2005*paciHost2005)*oecdABC +
                 lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                 ljiHost + pcHost + demHost +
                 lgdpHome + gdpgrowHome +
                 ljiHome + 
                 bit + comlang + colony + dist +
                 lassets + lage + lhosts +
                 (1|isocodeHome) + (1|isocodeHost) + (1 | naics2),
               data = orbis %>% filter(isocodeHost != "CN" & isocodeHome != "CN"),
               family = binomial(link="logit"), nAGQ = 0
)
summary(no.cn)
toc() # 28 seconds

tic()
no.count <- glmer(subsid ~ oecdABC * paciHost2005 +
                    I(paciHost2005*paciHost2005)*oecdABC +
                    lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                    ljiHost + pcHost + demHost +
                    lgdpHome + gdpgrowHome +
                    ljiHome + 
                    bit + comlang + colony + dist +
                    lassets + lage + lhosts +
                    (1|isocodeHome) + (1|isocodeHost) + (1 | naics2),
                  data = orbis %>% filter(isocodeHost != "IL" & isocodeHome != "IL" &
                                            isocodeHost != "ZA" & isocodeHome != "ZA"),
                  family = binomial(link="logit"), nAGQ = 0
)
summary(no.count)
toc() # 33 seconds

### Table E.2 ----
cf_map <- c("oecdABCI(paciHost2005 * paciHost2005)" = "OECD Ratifier X Host PACI sq.", 
            "oecdABCpaciHost2005" = "OECD Ratifier X Host PACI",
            "oecdABC" = "OECD Ratifier",
            "I(paciHost2005 * paciHost2005)" = "Host PACI sq.",
            "paciHost2005" = "Host PACI",
            "ljiHost" = "Host Judiciary Indep.", 
            "ljiHostljiHome" = "Host Judiciary Indep. X Home Judiciary Indep.",
            "ljiHome" = "Home Judiciary Indep.",
            "lgdpHost" = "Host GDP (log)", 
            "gdppcHost" = "Host GDP per capita",
            "fdigdpHost" = "Host FDI (GDP perc.)", 
            "tradeHost" = "Host Trade (GDP perc.)", 
            "pcHost" = "Host POLCON III", 
            "demHost" = "Host Democracy",
            "lgdpHome" = "Home GDP (log)", 
            "gdpgrowHome" = "Home GDP Growth (perc.)", 
            "bit" = "Dyad BIT", 
            "comlang" = "Dyad Common Language", 
            "colony" = "Dyad Colonial Relation", 
            "dist" = "Dyad Distance",
            "lassets" = "Firm Assets (log)", 
            "lage" = "Firm Age (log)", 
            "lhosts" = "Firm Host Countries (log)",
            "(Intercept)" = "Constant")

rows <- tribble(~term, ~original.BB, ~no.cn, ~no.count, 
                "Random intercepts", "Yes", "Yes", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(original.BB@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.cn@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.count@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(original.BB@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.cn@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.count@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(original.BB)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.cn)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.count)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(49:51, 53)

modelsummary(list("(1)" = original.BB,
                  "(2)" = no.cn,
                  "(3)" = no.count),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/supp_tab_e2.tex",
             title = "Firm-level data. Further robustness tests of multilevel logit models"
             )

## Extended dv ----
tic()
base.ex <- glmer(subsidExtend ~ oecdABC * paciHost2005 +
                   oecdABC * I(paciHost2005*paciHost2005) +
                   (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
                 data = orbis,
                 family = binomial(link = "logit"), nAGQ = 0)
summary(base.ex)
toc() # 19 seconds

tic()
country.ex <- glmer(subsidExtend ~ oecdABC * paciHost2005 +
                      oecdABC * I(paciHost2005*paciHost2005) +
                      lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                      ljiHost + pcHost + demHost +
                      lgdpHome + gdpgrowHome +
                      ljiHome + 
                      (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
                    data = orbis,
                    family = binomial(link = "logit"), nAGQ = 0)
summary(country.ex)
toc() # 18 seconds

tic()
dyad.ex <- glmer(subsidExtend ~ oecdABC * paciHost2005 +
                   oecdABC * I(paciHost2005*paciHost2005) +
                   lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                   ljiHost + pcHost + demHost +
                   lgdpHome + gdpgrowHome +
                   ljiHome + 
                   bit + comlang + colony + dist +
                   (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
                 data = orbis,
                 family = binomial(link="logit"), nAGQ = 0)
summary(dyad.ex)
toc() # 28 seconds

tic()
full.ex <- glmer(subsidExtend ~ oecdABC * paciHost2005 +
                   oecdABC * I(paciHost2005*paciHost2005) +
                   lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                   ljiHost + pcHost + demHost +
                   lgdpHome + gdpgrowHome +
                   ljiHome + 
                   bit + comlang + colony + dist +
                   lassets + lage + lhosts +
                   (1|isocodeHome) + (1|isocodeHost) + (1|dyad),
                 data = orbis,
                 family = binomial(link="logit"), nAGQ = 0)
summary(full.ex)
toc() # 33 seconds

tic()
full.ind.ex <- glmer(subsidExtend ~ oecdABC * paciHost2005 +
                       oecdABC * I(paciHost2005*paciHost2005) +
                       lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                       ljiHost + pcHost + demHost +
                       lgdpHome + gdpgrowHome +
                       ljiHome + 
                       bit + comlang + colony + dist +
                       lassets + lage + lhosts +
                       (1|isocodeHome) + (1|isocodeHost) + (1|naics2)  + (1|dyad),
                     data = orbis,
                     family = binomial(link="logit"), nAGQ = 0)
summary(full.ind.ex)
toc() # 41 seconds

### Table E.3 ----
cf_map <- c("oecdABCI(paciHost2005 * paciHost2005)" = "OECD Ratifier X Host PACI sq.", 
            "oecdABCpaciHost2005" = "OECD Ratifier X Host PACI",
            "oecdABC" = "OECD Ratifier",
            "I(paciHost2005 * paciHost2005)" = "Host PACI sq.",
            "paciHost2005" = "Host PACI",
            "lgdpHost" = "Host GDP (log)", 
            "gdppcHost" = "Host GDP per capita",
            "fdigdpHost" = "Host FDI (GDP perc.)", 
            "tradeHost" = "Host Trade (GDP perc.)", 
            "ljiHost" = "Host Judiciary Indep.", 
            "pcHost" = "Host POLCON III", 
            "demHost" = "Host Democracy",
            "lgdpHome" = "Home GDP (log)", 
            "gdpgrowHome" = "Home GDP Growth (perc.)", 
            "ljiHome" = "Home Judiciary Indep.",
            "bit" = "Dyad BIT", 
            "comlang" = "Dyad Common Language", 
            "colony" = "Dyad Colonial Relation", 
            "dist" = "Dyad Distance",
            "lassets" = "Firm Assets (log)", 
            "lage" = "Firm Age (log)", 
            "lhosts" = "Firm Host Countries (log)",
            "(Intercept)" = "Constant")

rows <- tribble(~term, ~base.ex, ~country.ex, ~dyad.ex, ~full.ex, ~full.ind.ex,
                "Random intercepts", "Yes", "Yes", "Yes", "Yes", "Yes",
                "Industry intercepts", "", "", "", "", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(base.ex@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(country.ex@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(dyad.ex@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.ex@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(full.ind.ex@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(base.ex@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(country.ex@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(dyad.ex@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.ex@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(full.ind.ex@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(base.ex)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(country.ex)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(dyad.ex)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.ex)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(full.ind.ex)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(47:50, 52)

modelsummary(list("(1)" = base.ex,
                  "(2)" = country.ex,
                  "(3)" = dyad.ex,
                  "(4)" = full.ex,
                  "(5)" = full.ind.ex),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/supp_tab_e3.tex",
             title = "Firm-level data. The effect of the OECD Convention on probability of subsidiary incorporation. Multilevel logit models. Extended data"
)

# Placebos ----
## Host country ----
tic()
no.bribe.base.host <- glmer(subsid ~ oecdABC * paciHost2005 +
                              oecdABC * I(paciHost2005*paciHost2005) +
                              (1|isocodeHome) + (1|isocodeHost) + (1|dyad) + (1 | naics2),
                            data = orbis %>% filter(host_bribery == 0),
                            family = binomial(link="logit"), nAGQ = 0
)
summary(no.bribe.base.host)
toc() # 7 seconds

tic()
no.bribe.full.host <- glmer(subsid ~ oecdABC * paciHost2005 +
                              oecdABC * I(paciHost2005*paciHost2005) +
                              lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                              ljiHost + pcHost + demHost +
                              lgdpHome + gdpgrowHome +
                              ljiHome + 
                              bit + comlang + colony + dist +
                              lassets + lage + lhosts +
                              (1|isocodeHome) + (1|isocodeHost) + (1 | naics2)  + (1|dyad),
                            data = orbis %>% filter(host_bribery == 0),
                            family = binomial(link="logit"), nAGQ = 0
)
summary(no.bribe.full.host)
toc() # 9 seconds

tic()
yes.bribe.base.host <- glmer(subsid ~ oecdABC * paciHost2005 +
                               oecdABC * I(paciHost2005*paciHost2005) +
                               (1|isocodeHome) + (1|isocodeHost) + (1|dyad) + (1 | naics2),
                             data = orbis %>% filter(host_bribery == 1),
                             family = binomial(link="logit"), nAGQ = 0
)
summary(yes.bribe.base.host)
toc() # 17 seconds

tic()
yes.bribe.full.host <- glmer(subsid ~ oecdABC * paciHost2005 +
                               oecdABC * I(paciHost2005*paciHost2005) +
                               lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                               ljiHost + pcHost + demHost +
                               lgdpHome + gdpgrowHome +
                               ljiHome + 
                               bit + comlang + colony + dist +
                               lassets + lage + lhosts +
                               (1|isocodeHome) + (1|isocodeHost) + (1 | naics2)  + (1|dyad),
                             data = orbis %>% filter(host_bribery == 1),
                             family = binomial(link="logit"), nAGQ = 0
)
summary(yes.bribe.full.host)
toc() # 40 seconds

### Table F.1 ----
rows <- tribble(~term, ~yes.bribe.base.host, ~yes.bribe.full.host, ~no.bribe.base.host, ~no.bribe.full.host,
                "Random intercepts", "Yes", "Yes", "Yes", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(yes.bribe.base.host@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(yes.bribe.full.host@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.bribe.base.host@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.bribe.full.host@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(yes.bribe.base.host@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(yes.bribe.full.host@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.bribe.base.host@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.bribe.full.host@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(yes.bribe.base.host)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(yes.bribe.full.host)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.bribe.base.host)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.bribe.full.host)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(47:49, 51)

modelsummary(list("Test (1)" = yes.bribe.base.host,
                  "Test (2)" = yes.bribe.full.host,
                  "Placebo (3)" = no.bribe.base.host,
                  "Placebo (4)" = no.bribe.full.host),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/supp_tab_f1.tex",
             title = "Firm-level data. Placebo test: Effects of the OECD Convention in host countries that experienced pre-2005 enforcement events"
)

## Home country ----
tic()
no.bribe.base.home <- glmer(subsid ~ oecdABC * paciHost2005 +
                              oecdABC * I(paciHost2005*paciHost2005) +
                              (1|isocodeHome) + (1|isocodeHost) + (1|dyad) + (1 | naics2),
                            data = orbis %>% filter(home_bribery == 0),
                            family = binomial(link="logit"), nAGQ = 0
)
summary(no.bribe.base.home)
toc() # 3 seconds

tic()
no.bribe.full.home <- glmer(subsid ~ oecdABC * paciHost2005 +
                              oecdABC * I(paciHost2005*paciHost2005) +
                              lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                              ljiHost + pcHost + demHost +
                              lgdpHome + gdpgrowHome +
                              ljiHome + 
                              bit + comlang + colony + dist +
                              lassets + lage + lhosts +
                              (1|isocodeHome) + (1|isocodeHost) + (1 | naics2)  + (1|dyad),
                            data = orbis %>% filter(home_bribery == 0),
                            family = binomial(link="logit"), nAGQ = 0
)
summary(no.bribe.full.home)
toc() # 7 seconds

tic()
yes.bribe.base.home <- glmer(subsid ~ oecdABC * paciHost2005 +
                               oecdABC * I(paciHost2005*paciHost2005) +
                               (1|isocodeHome) + (1|isocodeHost) + (1|dyad) + (1 | naics2),
                             data = orbis %>% filter(home_bribery == 1),
                             family = binomial(link="logit"), nAGQ = 0
)
summary(yes.bribe.base.home)
toc() # 26 seconds

tic()
yes.bribe.full.home <- glmer(subsid ~ oecdABC * paciHost2005 +
                               oecdABC * I(paciHost2005*paciHost2005) +
                               lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                               ljiHost + pcHost + demHost +
                               lgdpHome + gdpgrowHome +
                               ljiHome + 
                               bit + comlang + colony + dist +
                               lassets + lage + lhosts +
                               (1|isocodeHome) + (1|isocodeHost) + (1 | naics2)  + (1|dyad),
                             data = orbis %>% filter(home_bribery == 1),
                             family = binomial(link="logit"), nAGQ = 0
)
summary(yes.bribe.full.home)
toc() # 40 seconds

### Table F.2 ----
rows <- tribble(~term, ~yes.bribe.base.home, ~yes.bribe.full.home, ~no.bribe.base.home, ~no.bribe.full.home,
                "Random intercepts", "Yes", "Yes", "Yes", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(yes.bribe.base.home@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(yes.bribe.full.home@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.bribe.base.home@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.bribe.full.home@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(yes.bribe.base.home@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(yes.bribe.full.home@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.bribe.base.home@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.bribe.full.home@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(yes.bribe.base.home)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(yes.bribe.full.home)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.bribe.base.home)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.bribe.full.home)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(47:49, 51)

modelsummary(list("Test (1)" = yes.bribe.base.home,
                  "Test (2)" = yes.bribe.full.home,
                  "Placebo (3)" = no.bribe.base.home,
                  "Placebo (4)" = no.bribe.full.home),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/supp_tab_f2.tex",
             title = "Firm-level data. Placebo test: Effects of the OECD Convention in home countries that experienced pre-2005 enforcement events"
)

## Industry ----
tic()
no.bribe.base <- glmer(subsid ~ oecdABC * paciHost2005 +
                         oecdABC * I(paciHost2005*paciHost2005) +
                         (1|isocodeHome) + (1|isocodeHost) + (1|dyad) + (1 | naics2),
                       data = orbis %>% filter(industry_bribery == 0),
                       family = binomial(link="logit"), nAGQ = 0
)
summary(no.bribe.base)
toc() # 8 seconds

tic()
no.bribe.full <- glmer(subsid ~ oecdABC * paciHost2005 +
                         oecdABC * I(paciHost2005*paciHost2005) +
                         lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                         ljiHost + pcHost + demHost +
                         lgdpHome + gdpgrowHome +
                         ljiHome + 
                         bit + comlang + colony + dist +
                         lassets + lage + lhosts +
                         (1|isocodeHome) + (1|isocodeHost) + (1 | naics2)  + (1|dyad),
                       data = orbis %>% filter(industry_bribery == 0),
                       family = binomial(link="logit"), nAGQ = 0
)
summary(no.bribe.full)
toc() # 9 seconds

tic()
yes.bribe.base <- glmer(subsid ~ oecdABC * paciHost2005 +
                          oecdABC * I(paciHost2005*paciHost2005) +
                          (1|isocodeHome) + (1|isocodeHost) + (1|dyad) + (1 | naics2),
                        data = orbis %>% filter(industry_bribery == 1),
                        family = binomial(link="logit"), nAGQ = 0
)
summary(yes.bribe.base)
toc() # 33 seconds

tic()
yes.bribe.full <- glmer(subsid ~ oecdABC * paciHost2005 +
                          oecdABC * I(paciHost2005*paciHost2005) +
                          lgdpHost + gdppcHost + fdigdpHost + tradeHost +
                          ljiHost + pcHost + demHost +
                          lgdpHome + gdpgrowHome +
                          ljiHome + 
                          bit + comlang + colony + dist +
                          lassets + lage + lhosts +
                          (1|isocodeHome) + (1|isocodeHost) + (1 | naics2)  + (1|dyad),
                        data = orbis %>% filter(industry_bribery == 1),
                        family = binomial(link="logit"), nAGQ = 0
)
summary(yes.bribe.full)
toc() # 27 seconds

### Table F.3 ----
rows <- tribble(~term, ~yes.bribe.base, ~yes.bribe.full, ~no.bribe.base, ~no.bribe.full,
                "Random intercepts", "Yes", "Yes", "Yes", "Yes",
                # save number of home and host countries to properly display it in the table
                "N. of host countries",  
                (length(unique(yes.bribe.base@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(yes.bribe.full@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.bribe.base@flist$isocodeHost)) +1) %>% as.character(),
                (length(unique(no.bribe.full@flist$isocodeHost)) +1) %>% as.character(),
                "N. of home countries",
                (length(unique(yes.bribe.base@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(yes.bribe.full@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.bribe.base@flist$isocodeHome)) +1) %>% as.character(),
                (length(unique(no.bribe.full@flist$isocodeHome)) +1) %>% as.character(),
                "Log Likelihood",
                summary(yes.bribe.base)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(yes.bribe.full)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.bribe.base)$log %>% as.numeric() %>% round(digits = 3) %>% as.character(),
                summary(no.bribe.full)$log %>% as.numeric() %>% round(digits = 3) %>% as.character())

attr(rows, "position") <- c(47:49, 51)

modelsummary(list("Test (1)" = yes.bribe.base,
                  "Test (2)" = yes.bribe.full,
                  "Placebo (3)" = no.bribe.base,
                  "Placebo (4)" = no.bribe.full),
             stars = TRUE,
             coef_map = cf_map,
             gof_map = c("nobs", "aic"),
             add_rows = rows,
             output = "tables/supp_tab_f3.tex",
             title = "Firm-level data. Placebo test: Effects of the OECD Convention in industries that experienced pre-2005 enforcement events"
)

#====# The End #====#